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Numerical  Simulation  of  Thunderstorm  Gust  Fronts 


1.  INTRODUCTION 

A  relatively  small  but  very  important  feature  of  thunderstorms  is  the  gust 
front.  The  gust  front  is  the  leading  edge  of  the  cold  downdraft  air  that  spreads 
out  upon  reaching  the  ground.  As  shown  in  Figure  1,  the  gust  front  acts  as  a 
wedge  that  lifts  low  level  air  into  the  updraft  of  the  thunderstorm.  This  main¬ 
tains  the  thunderstorm  by  either  supplying  a  continual  inflow  of  potentially  un¬ 
stable  air  in  the  case  of  a  supercell  storm  or  initiating  new  updraft  cells  in  the 
case  of  a  multicell  storm.  *  Gust  fronts  that  have  propagated  well  away  from 
their  parent  thunderstorms  can  be  responsible  for  the  initiation  of  new  thunder¬ 
storm  activity,  especially  at  the  intersection  of  the  gust  front  and  some  other 
air  mass  boundary  or  front.  Recent  studies  have  indicated  that  the  gust  front 
may  also  play  a  role  in  tornado  genesis.  ^ 

(Received  for  publication  29  Nov  1983) 
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Figure  1.  Cross-section  Through  a  Supercell  Type  Thunderstorm 


Observational  studies  have  shown  that  the  cold  air  outflow  from  the  thunder - 
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storm  can  be  treated  as  a  density  current.  A  schematic  of  an  atmospheric 
density  current  is  shown  in  Figure  2.  The  current  is  characterized  by  an  ele¬ 
vated  head  region,  an  overhanging  nose,  and  a  turbulent  wake  at  the  rear  of  the 
head.  The  arc  cloud  shown  in  Figure  2  is  a  commonly  observed  feature  that  is 
the  result  of  moist  boundary  layer  air  being  lifted  above  its  level  of  condensation 
during  the  passage  of  the  current. 

5.  Fritsch,  J.  M.  (1975)  Cumulus  dynamics:  Local  compensating  subsidence 
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Figure  2.  Schematic  Representation  of  Atmospheric 
Density  Current* 
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There  have  been  many  laboratory  studies  of  density  current  motion  and 

several  numerical  studies.  While  one  of  the  laboratory  studies  dealt  with 

g 

the  effect  of  ambient  wind  on  the  gust  front  motion,  none  of  these  studies  in¬ 
cluded  the  effect  of  the  environmental  shear  commonly  associated  with  severe 
storms  or  moist  processes  such  as  arc  cloud  formation.  These  studies  also 
did  not  address  the  question  of  "gust  front  decoupling",  the  process  by  which 

the  gust  front  separates  from  the  parent  storm  and  leads  to  the  storms  rapid 

.  8 
decay. 

This  research  is  meant  to  extend  the  work  of  the  previous  numerical  and 
laboratory  studies  through  the  use  of  a  numerical  model.  Specifically,  the 
objectives  were: 

(a)  Determine  the  effect  of  the  ambient  winds  on  gust  front  motion; 

(b)  Examine  the  effect  of  arc  cloud  generation  on  the  gust  front  speed  and 
intensity; 

(c)  Investigate  the  decoupling  phenomena  associated  with  gust  fronts; 

(d)  Determine  the  profile  of  lifting  that  is  associated  with  gust  frontal 
passage; 

9.  Simpson,  J.E.,  and  Britter,  R.E.  (1980)  A  laboratory  model  of  an  atmos¬ 
pheric  mesofront,  Quart.  J.  Roy.  Meteorol.  Soc,  106:485-300. 

10.  Daly,  B.J.,  and  Pracht,  W.E.  (1968)  Numerical  study  of  density  current 

surges,  Phys.  Fluids  10:15-30. 

11.  Mitchell,  K.E.,  and  Hovermale,  J.  B.  (1977)  A  numerical  investigation  ol 

the  severe  thunderstorm  gust  front,  Mon.  Wea.  Rev.  105:657-673. 

12.  Thorpe,  A.J. ,  Miller,  M.J. ,  and  Moncrieff,  M.W.  (1980)  Dynamical  tv  o- 

dimensional  downdraughts,  Quart.  J.  Roy.  Meteorol.  Soc.  106:463-48-4. 


(e)  Show  that  the  gust  front  can  be  accurately  represented  in  a  numerical 
model  that  uses  the  same  resolution  as  current  thunderstorm  models. 

One  aspect  of  the  study  pertaining  to  objectives  (a),  (b),  and  (c)  was  the 
development  of  a  more  general  gust-front  speed  equation  than  has  been  used  in 
previous  studies.  This  important  result  will  be  discussed  in  the  following  section. 


2.  THE  GUST  FRONT  SPFF.U  EQUATION 

A  major  concern  of  many  observational  and  laboratory  studies  of  gust  front 

and  density  current  motion  has  been  the  determination  of  an  equation  relating 
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the  speed  of  the  front  to  other  measurable  quantities.  *■  Most  of  these  studies 

13 

used  an  equation  of  the  form 

eh  .  (1) 

In  I’(|.  11',  V  is  th<  velocity  of  the  front,  g  the  acceleration  of  gravity,  h  the 
height  of  th*  denser  fluid  outflow  ,  Ae  the  difference  in  density  between  the  cur¬ 
rent  and  the  environment,  and  the  density  of  the  environment.  The  constant, 
k  ,  is  tin  internal  Inmile  number.  Htnjum  in^ showed  that  for  two  immiscible, 
mv  is.  Hi  fluids ,  k  -  \  2  . 

lii  1  q.  (li,  tin  hi  :ghl,  h,  refers  to  tin  depth  of  the  fluid  behind  the  head. 
Mmi  tin  ratio  of  the  h«  ight  ol  the  head  and  tfie  fluid  depth  behind  the  head  is 
nearly  *  onstant,  '  i  rt  definition  ol  th*  1  i  oud*  numln  r  allows  the  use  of  the 
height  of  th*  head,  H.  If  a*  assume  that  th*  tensity  ol  an  varies  with  tempera  - 
tui a  *  ,  -  [1  *  II  -  I  l  I  |  ,  a  hi  i  •  I  is  tin  vutuul  temperature  and,,  is 

u  Vo  V  v  O 

th*  dcnsitv  at  1  ,  then  Id-  ill  mav  n*  i<  *.  i  itt*  i.  ho  atmosphere  density  iur- 

rt  nt.s  a  .- 


A  1 

\  k  gll  — 

In  Eq.  (2),  AT  is  the  differ*  te  «  ;t.  v;i  to.*.  ’«  ••  p»  i  >t..i  *  *n,  .*.  , i  •• 

and  cold  air,  T  is  the  virtual  temper-. tu : «  *:  *h.  ■  ■  :  :  .  o  '  •  .-  . . - 

fined  Froude  number  fk  =  ^  h  it  k  ». 


CM 


13.  '  -Mamin,  J.  B.  119681  firavit. 

•ch.  31:209-248. 
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Although  theoretically  the  value  of  k*  =  \.  2  and  k  «  1.0,  the  measured  values 

of  these  quantities  in  the  atmosphere  and  laboratory  are  substantially  smaller. 

This  is  because  the  actual  flows  occur  in  fluids  that  have  viscous  and  surface 

drag  effects.  These  effects  tend  to  reduce  the  speed  from  the  theoretical  limit, 

and  this  is  incorporated  into  the  value  of  the  Froude  number.  There  has  been  a 

fairly  wide  interval  of  k  values  obtained  for  atmospheric  density  currents  and 

8 

gust  fronts  ranging  from  k  =  0.  7  0  to  k  =  1.  08  . 

The  driving  force  of  the  density  current  is  the  pressure  gradient  force 
resulting  from  the  increased  hydrostatic  pressure  in  the  cold  air.  This  is 
reflected  in  Eqs.  (1)  and  (2),  which  may  be  rewritten  in  terms  of  the  pressure. 

Since  the  difference  in  surface  hydrostatic  pressure  between  the  head  region 
and  the  environment  is  Ap  =  gHAp,  Eq.  (1)  may  be  written  (again  with  the  rede¬ 
fined  Froude  number) 

v-[*r  • 

where  Ap  is  the  difference  in  surface  hydrostatic  pressure  between  the  density 
current  head  and  the  environment. 

For  the  constant  density  fluids  assumed  in  the  derivation  of  Eq.  (1),  all 
three  of  these  equations  will  give  identical  results.  However,  for  atmospheric 
density  currents,  the  results  given  by  these  three  equations  will  differ  substan¬ 
tially.  This  is  because  Eqs.  (1)  and  (2)  deal  only  with  the  density  differences 
over  the  height  of  the  outflow,  while  Eq.  (3)  includes  the  integrated  effect  of  the 
entire  depth. 

In  applying  Eq.  (3)  to  observed  atmospheric  density  currents,  some  assump¬ 
tions  must  be  made.  The  pressure  difference,  Ap,  refers  to  the  surface  hydro¬ 
static  pressure  in  Eq.  (3),  but  the  measured  pressure  is  the  total  surface  pres- 
8 

sure.  Wakimoto  has  shown,  however,  that  the  maximum  surface  pressure 
increase  measured  during  the  passage  of  the  head  is  primarily  hydrostatic  and 
can  therefore  be  used  directly  in  Eq.  (3).  The  numerical  results  of  Mitchell  and 
Hovermale11  also  indicate  this.  The  density,  p,  in  Eq.  (3)  is  the  average  density 
of  the  warm  air  over  the  depth  of  the  current,  H.  For  ease  of  application,  how¬ 
ever,  the  surface  density  may  be  used  in  Eq.  (3)  with  virtually  no  error. 

It  is  common  for  gust  fronts  to  propagate  through  environments  that  have  an 

ambient  wind.  The  component  of  the  wind  parallel  to  the  gust  front  motion  will 

o 

affect  its  speed  according  to 


v  ■  kbN 


+  0.  G2U  , 


(4) 


where  U  is  the  wind  parallel  to  the  gust  front  motion  averaged  over  the  depth  of 
the  head  and  positive  in  the  direction  of  gust  front  motion. 

In  order  to  determine  the  value  of  k  that  applies  to  gust  fronts  when  the 
pressure  form  of  the  equation  is  used,  Eqs.  (3)  and  (4)  were  applied  to  2  0  ob¬ 
served  gust  fronts.  These  included  the  case  studied  by  Charba,  r>  the  three  cases 
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of  Wakimoto,  and  the  type  I,  II,  and  III  cases  of  Goff.  The  data  on  these  cases 
included  the  recorded  surface  pressure  rise,  gust  front  speed,  and  ambient 
winds  information.  The  surface  value  of  the  density  in  the  warm  air  was  used 
in  all  cases. 

Figure  3  shows  these  observations  plotted  so  that  they  should  fall  on  a 
straight  line  with  k  as  the  slope.  The  best  fit  line  through  the  origin  gives 
k  =  0.79.  The  correlation  coefficient  is  0.  84.  The  numerical  results  of 
Mitchell  and  Hovermale1*  are  also  plotted  in  Figure  3,  but  were  not  used  in 
determining  the  best  fit  line. 
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Figure  3.  Observations  of  Gust  Front  Speed  vs  Surface 
Pressure  Rise 
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The  numerical  results  shown  in  Figure  3  fit  the  observational  results  very 
well.  This  shows  that  the  numerical  approach  taken  by  Mitchell  and  Hovermale 
captures  the  essential  features  of  atmospheric  density  current  motion,  and 
justifies  the  extension  of  this  approach  taken  in  the  present  study. 

3.  THE  NUMERICAL  MODEL 


3. 1  The  Equations 

The  model  used  to  obtain  the  numerical  results  reported  in  this  study  is  an 
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improved  version  of  the  2-D,  Boussinesq,  moist  model  developed  by  Seitter. 
The  domain  is  40  km  in  the  horizontal  and  10  km  in  the  vertical  in  the  x-z  plane 
with  a  grid  spacing  of  500  m  in  both  directions.  The  equations  governing  the 
vorticity,  r),  temperature,  T,  water  vapor  mixing  ratio,  q,  and  liquid  water 
mixing  ratio,  m,  are,  respectively, 


If  +  u|2-+  w|a  - 
9t  9x  9z 


1  9Tv  „  9m  .  K  2 

e:  "55T  +  g  ar  +  Kv  " 


9T  ,  9T  9T  ,  _  L  „  .  -.2™ 

inr  +  u9T  +  w si  +  wrd  =  -G  +  Kv  T  * 

p 


!f  +  u!f +  wH-  =  -G  +  > 


9m  ,  9m  ,  9m  .  9  _ 

_+u_„+w_+_  (vtm)  =  G  +  KV  m 


The  velocities,  u  and  w,  are  written  in  terms  of  a  streamfunction,  i  , 


-  h. 

'  9z  ’ 


9u  9w  2 

*  =  97  '  aT  =  v  1 


In  the  above. 


14.  Seitter,  K.  L.  (1082)  The  Dynamical  Structure  of  Squall-line  Type  Thunder¬ 
storms,  unpublished  Ph.  D.  thesis.  Dept.  Gcophvs.  Sci.  ,  Univ.  of 
Chicago,  130  pp. 


T  *  (1  +  0.  61q)T  =  virtual  temperature, 
eo  =  300K, 

I'd  =  9.76°C/km  =•  dry  adiabatic  lapse  rate, 

L  =  latent  heat  of  condensation, 

Cp  *  specific  heat  at  constant  pressure  for  dry  air, 

G  =  liquid  water  generation  term, 

V.  =  precipitation  fallspeed, 

K  =  eddy  viscosity. 

The  equations  are  finite  differenced  using  the  Arakawa  Jacobian  for  the 

advective  terms  and  the  consistent  nine-point  finite  difference  form  for  the  non- 
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advective  terms.  Time  integration  is  performed  using  the  Adams-Bashforth 

method  with  diffusion.  *5  The  Poisson  equation  for  the  streamfunction,  Eq.  (10), 

is  solved  after  each  time  step  by  successive -over-relaxation.  This  equation  is 

finite  differenced  using  a  fourth  order  differencing  form  to  improve  accuracy. 

The  condensation  scheme  that  determines  the  value  of  G  at  each  grid  point 

is  described  in  Seitter  and  Kuo.  15  Water  condensed  up  to  lg/kg  is  assumed  to 

be  cloud  water,  and  beyond  this  value  it  is  assumed  to  fall  as  precipitation.  The 

1 

variable  fallspeed  formulation  follows  that  of  Liu  and  Orville, 

(  0  for  m  £  1.  0 

vt-  <  ,  (11) 

ls.32(m  -  1.0)0,2  form  >1.0 
for  m  in  g/kg  and  V  in  m/sec. 

3.2  Surface  Drag  and  Eddy  Diffusion  Formulation 

The  surface  is  treated  as  a  "limited  slip"  surface  following  Mitchell  and 
Hovermale.  **  In  this  formulation  the  boundary  condition  is  applied  as  freeslip 
with  a  surface  drag.  Mitchell  and  Hovermale**  applied  this  condition  in  a 
primitive  equation  model.  Its  application  in  a  vortieity  model  is  less  obvious 
and  requires  careful  consideration  for  consistency.  The  surface  drag  is  applied 
as  though  it  retarded  only  the  surface  velocity,  which  in  turn  modifies  the 

15.  Seitter,  K.  L.  ,  and  Kuo,  H.  -L.  (1983)  The  dynamical  structure  of  squall¬ 

line  type  thunderstorms,  J.  Atmos.  Sci.  40:283 1  -2854. 

16.  Liu,  J.Y.,  and  Orville,  H.D.  (1969)  Numerical  modeling  of  precipitation  and 

cloud  shadow  effects  on  mountain  cumuli,  J.  Atmos.  Sci.  26:1282-1298. 
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vorticity  at  the  grid  level  just  above  the  surface.  The  surface  vorticity  is  held 
at  zero  as  implied  by  a  free-slip  condition. 

To  evaluate  the  surface  drag  effect  we  must  first  calculate  the  surface 
velocity  based  on  the  interior  flow.  Using  the  definition  of  the  vorticity  and 
Eq.  (P)  we  may  write 


1 


(12) 


Then,  using  a  centered  difference  for  9u/9z,  we  can  write  the  velocity  at  the 
surface  (level  1),  at  the  horizontal  grid  point  i,  in  terms  of  the  interior  quantities 


u 


i,  1 


i.  3 


(2Az) 


(13) 


In  general,  the  surface  velocity  calculated  in  this  way  will  yield  u.  „  -  u. 

1*  “  1#  A 

f  0,  that  is,  Au/Az  }  0.  This  is,  however,  not  a  violation  of  the  free-slip  con¬ 
dition.  The  analytic  free-slip  condition  is  9u/9z  =  0  at  the  surface.  This  leads 
to  the  specification  of  n  =  0  on  the  lower  boundary,  since  w  =  9w/9x  =  0  there. 
The  specification  of  n  =  0  at  the  surface  satisfies  the  analytic  free-slip  condition 
without  reference  to  the  velocity  at  level  2. 

We  now  allow  the  surface  velocity  calculated  in  Eq.  (13)  to  be  retarded  by 

a  surface  drag.  The  deceleration  of  u.  .  by  a  bulk  aerodynamic  drag  may  be 
11  1,1 

written 


3u. 


i.  1 


~W 


drag 


c  .  u.  . 

d  i.l 
Az 


LA.  - 

1.  1 


(14) 


where  c^  is  the  drag  coefficient  and  Az  represents  the  depth  of  the  layer  over 
which  the  drag  acts  (one  vertical  grid  spacing  in  this  case). 

The  surface  vorticity  is  unaffected  by  the  deceleration  of  u.  .,  but  the 
vorticity  at  level  2  must  be  modified.  The  rate  of  change  of  the  vorticity  at 
level  2  due  to  the  drag  is  clearly 


an, 


ij  2 


at 


drag 


9 

9t 


r  -in.  i 

l  9z  dragj 


(la) 


Since  the  drag  is  only  applied  at  level  1,  Eq.  (15)  reduces  to 


The  surface  drag  effects  are  then  included  in  the  model  by  replacing  the 

2 

KV  r)  term  in  Eq.  (5)  by  the  right  side  of  Eq.  (16)  in  the  calculation  of  the  time 
rate  of  change  of  vorticity  at  level  2.  The  diffusion  term  is  omitted  at  this  level 
because  it  would  tend  to  diffuse  zero  vorticity  from  the  surface  into  the  interior 
of  the  domain,  which  is  not  desirable.  If  we  set  c^  =  0,  we  obtain  a  true  free- 
slip  boundary  condition.  In  the  model  simulations,  the  drag  coefficient  was  set 
at  c^  =  0.  02  following  Mitchell  and  Hovermale.  ^ 

Note  that  the  surface  velocity  calculated  in  Eq.  (13)  is  a  "true"  surface 
velocity  and  not  a  nonphysical  numerical  quantity.  Therefore,  this  velocity  can 
be  used  to  advect  surface  values  of  the  trace  variables  with  physical  justification. 
This  formulation  can  easily  be  generalized  to  an  anelastic  model. 

The  diffusion  of  quantities  by  turbulent  eddies  is  treated  with  a  Fickian 
diffusion  term  using  a  nonlinear  eddy  viscosity,  K.  We  let 


K  =  Kq  +  K*  |  v  n  I 


where  K  and  K*  are  constants.  K  is  a  constant  eddy  coefficient  that  is  re- 
quired  to  stabilize  the  Adams-Bashforth  method.  A  value  of  at  least  150  m  /sec 
is  required  to  allow  a  relatively  large  timestep  (10  sec),  so  that  the  model  can 
be  economically  used  for  a  large  number  of  simulations.  Two  forms  of  K  were 
used  in  the  simulations  presented  here.  These  will  be  discussed  in  connection 
with  the  two  cold  air  source  types  presented  later  in  this  section. 

3.3  Other  Boundary  Conditions 

Lateral  boundary  conditions  are  always  an  important  concern  in  limited 
area  modeling.  This  model  uses  a  2-D  radiation  condition  on  the  vorticity 


3m  9>v  9  q, 

_ b  .  r  _ b  .  r  _ b.  _  n 

at  x  ax  z  az  ~  u 


where  is  the  value  of  the  vorticity  at  the  boundary.  In  Eq.  (18),  and  Cz 
represent  the  speed  of  an  outgoing  wavelike  disturbance  perpendicular  to  and 
parallel  to  the  boundary,  respectively,  and  are  calculated  from  the  interior 
solution.  The  numerical  implementation  of  this  condition  is  outlined  in  Seitter 
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and  Kuo.  The  other  variables  satisfy 


u  into  domain 


(19a 


i£  =  o 

at 

+  u  =0  u  out  of  domain  ,  (19b! 

where  <j>  =  T,q,m,  and  u  =  dtl /dz  at  the  lateral  boundary.  Neumann  conditions, 
a,/  /3x  =  0,  are  specified  for  the  streamfunction. 

The  top  and  bottom  boundaries  are  taken  as  rigid  plates,  so  the  stream  - 
function  is  constant  on  these  boundaries.  The  vorticity  is  set  equal  to  zero  on 
the  top  boundary.  The  other  variables  are  held  fixed  on  the  top  boundary  with 
the  temperature  and  water  vapor  mixing  ratio  held  equal  to  the  environmental 
values  and  the  liquid  water  specified  as  zero.  On  the  bottom  boundary,  rain  is 
permitted  to  fall  "through"  the  surface  by  specifying 


3m 

at 


+  ir  (vtm)  = 


(2o; 


The  bottom  is  allowed  to  cool  by  evaporation  if  rain  is  present  by  specifying 

i 


* 

t 

I 


U 


3T 
o  3x 


(2 1; 


!?■ +  uo  = 


-  G 


(22! 


where 


C  -E  m  >  0,  RH  <  90% 

1^0  otherwise  ,  (23) 


—  R  —  1 

where  E  =4.0/  10  sec  is  the  constant  evaporation  rate  used  in  the  model. 
Equations  (21)  and  (22)  allow  rain  to  evaporate  up  to  a  relative  humidity  of 
90  percent,  and  for  the  surface  values  of  T  and  q  to  be  advected  by  the  surface 
wind.  The  surface  wind,  UQ,  is  given  by  Eq.  (13). 


3.4  Cold  Air  Source 

The  density  current  is  initiated  by  a  cold  air  source  located  near  the  left 
boundary.  Two  types  of  cold  air  sources  were  used.  The  first,  referred  to 
as  M-type,  is  a  fixed  temperature  source  based  on  the  work  of  Mitchell  and 
Hovermale.  * 1  For  this  source,  the  temperature  in  the  source  region  is  given 
initially  by 

T  =  Tenv.  -  ^T^°s<2ir/L  (z-zQ))  cos(2t/l  (x-Xq))]  ,  (24) 

for  |z  -  Zq|  s.  L/4  and  |x  -  Xq|  ^  L/4,  where  L  =  12  km.  Tenv  is  the  temper¬ 
ature  of  the  undisturbed  environment  at  that  height.  Once  the  integration  is 
begun  only  the  top  half  of  this  region  is  maintained  at  this  temperature.  While 
this  source  has  the  desirable  property  of  providing  cold  air  at  a  fixed  tempera- 
ture,  it  requires  a  large  eddy  viscosity,  K  =  000  m  /sec,  to  suppress  noise 

that  is  generated  near  the  source.  In  the  eddy  viscosity  formulation,  K  = 

2  ° 

900  m  /sec  and  K*  =  0,  so  for  the  simulations  using  the  M-type  source,  the 

formulation  reduces  to  a  constant  eddy  viscosity.  The  model  was  stable  for 

smaller  values  of  K,  but  the  value  used  confined  the  noise  to  the  immediate 

region  of  the  cold  air  source. 

A  second  type  of  cold  air  generator  used  is  based  on  the  work  of  Thorpe 
12 

et  al  and  is  referred  to  as  the  Q-type  source.  This  source  applies  a  constant 
cooling  rate  of  the  form 

4^-  =  -Q[cos(2tt/L  (z  -  z0) )  cos(2jr/L  (x -xq)  )]  ,  (25) 

for  |z  -  Zq|  S  L/4  and  |x  -  x^|  s.  L/4.  Q  is  the  cooling  rate  in  the  generator 
region.  This  source  produces  a  smoother  temperature  field  in  the  source 
region  and  allows  the  use  of  the  full  nonlinear  eddy  viscosity  formulation.  In 

the  Q-type  simulations,  K  =  150  m  /sec  and  K*  =  1.2  /  10ll)  m  .  With  these 

^  2 
values  of  K  and  K*,  K  reached  maximum  values  of  about  fiOO  m  /sec  in  the 

°  2 
turbulent  head  region  while  the  domain  average  never  exceeded  2  00  m  / sec. 

The  density  currents  produced  using  the  Q-type  source  in  the  present  study 

were  not  as  steady  as  those  obtained  with  the  M-type  source  because  the  flux 

of  cold  air  tended  to  pulsate. 
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4.  M  TYPE  SOURCE  RESULTS 


4.1  Stationary  Source 

This  section  will  describe  the  results  obtained  with  the  model  using  the 
M-type  cold  air  source  and  a  dry  atmosphere  with  no  environmental  winds.  In 
all  simulations,  the  atmosphere  had  a  constant  lapse  rate  and  a  surface  temper¬ 
ature  of  295  K.  The  results  for  the  M-type  simulations  are  summarized  in 
Table  1.  The  quantities  listed  in  this  table  were  measured  in  the  numerical 
results  after  the  density  current  became  steady.  Steadiness  was  generally 
reached  by  the  time  the  leading  edge  of  the  current  was  10  km  from  the  source 
region,  which  usually  occurred  by  15  minutes  into  the  integration. 

Figure  4a  shows  the  temperature  distribution  for  one  of  these  simulations, 
case  MD2,  after  20  min  of  integration.  Note  the  accurate  representation  of  the 
head  region  compared  to  Figure  2.  No  nose  is  evident  in  this  figure,  but 
observational  and  laboratory  studies  indicate  that  the  nose  would  be  below 
500  m  and  thus  not  resolvable  in  the  model.  The  outflow  depth  is  quite  large 

in  this  simulation  because  of  the  small  static  stability  Ofl/9z  =  0.  2°C/km),  but 
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is  almost  exactly  that  predicted  by  Mitchell  for  this  stability. 

The  streamfunction  for  case  MD2  is  analyzed  in  Figure  4b.  The  symmetry 
of  the  streamlines  around  the  source  region  shows  the  ability  of  the  2-D  lateral 
boundary  condition  to  accurately  allow  inflow  and  outflow  through  the  boundary. 

A  packing  of  the  streamlines  is  evident  at  the  surface  near  the  frontal  boundary 

of  the  outflow.  The  high  wind  gust  just  behind  the  front  associated  with  this 

6  7 

packing  is  a  commonly  observed  feature  during  gust  frontal  passage  '  and  has 
been  shown  to  be  a  result  of  the  surface  drag.  11 

A  more  realistic  depth  is  obtained  for  the  outflow  if  a  larger  and  more 
realistic  environmental  stability  is  used.  The  result  of  one  of  these  simula¬ 
tions,  case  MD9,  is  shown  in  Figure  5.  In  this  figure,  the  potential  tempera¬ 
ture  is  analyzed  instead  of  the  temperature.  In  this  representation,  the  outflow 
is  delineated  by  the  295  K  isentrope. 

4.2  Sensitivity  to  Model  Parameters 

The  eddy  viscosity  and  surface  drag  coefficient  are  two  parameters  that 
are  difficult  to  choose  on  purely  physical  grounds.  Therefore,  three  sensitivity 
simulations  were  carried  out  which  were  identical  to  MD2  except  for  changes 
in  one  of  these  parameters.  The  results  of  these  simulations  are  given  in 
Table  2. 

17.  Mitchell,  K.E.  (1975)  A  Numerical  Investigation  of  Severe  Thunderstorm 

Gust  Fronts,  unpublished  M.S.  thesis,  the  Pennsylvania  State  University, 
104  pp.  TaIso  available  from  NTIS,  No.  N70-17737  ^ 7GI ) 
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Table  2.  Sensitivity  Experiment  Results 


2 

In  case  MD2A,  the  eddy  viscosity  was  halved  to  a  value  of  K  =  450  m  /sec. 
The  resulting  density  current  was  nearly  identical  to  MD2,  but  propagated  a 
little  over  1  m/sec  faster  due  to  the  decreased  viscous  drag.  This  resulted  in 
an  increase  in  the  Froude  number  of  about  8  percent.  A  detailed  analysis  of  the 
effect  of  viscosity  on  the  speed  of  the  gust  front  was  carried  out  by  Daly  and 
Pracht.  10  Although  exact  comparison  with  their  results  is  not  possible,  they 
are  consistent  with  the  results  found  here. 

For  cases  MD2B  and  MD2C,  the  surface  drag  was  set  at  c^  =  0  and  = 

0.  04,  respectively,  instead  of  the  cd  =  0.  02  used  in  the  rest  of  the  simulations. 

As  would  be  expected,  the  gust  front  speed  was  inversely  related  to  the  surface 
drag  while  the  current  shape  and  surface  hydrostatic  pressure  rise  were  relative¬ 
ly  unaffected.  The  percentage  change  in  the  frontal  speed  for  these  simulations 
was  almost  identical  to  the  corresponding  experiments  in  the  more  complete 
sensitivity  tests  performed  by  Mitchell  and  Hovermale.  ** 

4.3  Moving  Source 

As  the  thunderstorm  propagates,  the  downdraft  represents  a  moving  cold 
air  source,  so  it  is  of  interest  to  consider  the  effect  that  movement  of  the  source 
has  on  the  resulting  density  current.  This  extension  to  previous  studies  is  made 
possible  by  the  excellent  lateral  boundary  conditions  used  in  the  current  model. 
Numerically,  it  is  much  easier  to  subtract  a  velocity  from  the  domain  than  to 
advect  the  source,  so  this  is  the  approach  taken.  Note  that  the  velocity  must 
also  be  subtracted  in  the  surface  drag  formulation  so  that  the  proper  drag  is 
applied. 

An  interesting  result  pertains  to  the  time  at  which  the  source  is  set  in 
motion.  If  the  density  current  is  allowed  to  propagate  on  the  order  of  15  km 
away  from  the  source  while  the  source  is  stationary,  and  then  the  source  is  set 
in  motion  at  a  speed  less  than  or  equal  to  the  dens  current  speed,  the  resulting 
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density  current  is  indistinguishable  from  one  in  which  the  source  is  stationary 
for  the  entire  integration.  In  fact,  if  the  source  is  set  in  motion  at  the  speed  of 
the  current,  the  source/current  system  remains  nearly  steady  with  the  current 
showing  virtually  no  propagation  relative  to  the  source.  In  this  situation  we  may 
say  that  the  gust  front  is  "decoupled"  from  the  source  in  the  sense  that  its  speed 
is  determined  by  its  local  characteristics. 

If  the  source  is  moving  at  the  start  of  the  integration,  the  resulting  density 
current  is  modified  relative  to  that  given  by  a  stationary  source.  This  is  illus¬ 
trated  by  comparing  the  results  in  Table  1  of  the  moving  source  cases  (MD3, 
MD4,  and  MD10)  with  the  corresponding  stationary  source  cases  (MD2  or  MDf). 
In  the  moving  source  simulations,  cold  air  is  "piled  up"  near  the  source  during 
the  early  stages  of  the  integration  before  the  current  propagates  a  sufficient 
distance  to  be  decoupled.  This  greater  depth  results  in  a  higher  hydrostatic 
pressure  and  thus  a  greater  speed. 

Figure  fi  shows  the  density  current  that  results  from  having  the  source 
move  at  10  m/sec.  This  simulation,  case  MD4,  was  identical  to  MD2  except 
for  the  movement  of  the  source.  Note  that  Figure  6  is  in  the  source  relative 
frame,  so  that  the  front  is  moving  faster  than  in  case  MD2  even  though  it  has 
traveled  less  far  in  the  domain  by  the  20-min  point.  The  piling  up  of  cold  air 
is  quite  evident  when  this  figure  is  compared  with  Figure  4a. 


Figure  6.  Temperature  Distribution  for  Case  MD4  at  20  min.  Source  rela¬ 
tive  winds  are  plotted  in  m/sec 


Since  most  severe  thunderstorms  develop  in  environments  exhibiting  vertical 
shear,  constant  vertical  shear  was  included  in  simulations  MDfi  and  MD7.  The 
results  of  case  MDG  after  20  min  of  integration  are  shown  in  Figure  7a.  For 

this  simulation  the  source  moved  at  10  m/sec  and  the  environmental  shear  was 

-3  - 1 

2/10  sec  .  Because  of  the  shear,  the  frontal  boundary  tends  to  be  retarded 


Figure  7a.  Temperature  Distribution  lor  Case  Ml)<>  at  20  min.  Source  relative 
winds  are  plotted  in  in  .'see 

more  near  the  surface  than  above  so  the  frontal  boundary  becomes  more  vertical 
The  current  shapes  in  these  simulations  were  very  similar  to  the  gust  front 

g 

observed  by  C'harba.  The  strcamfunetion  for  this  case,  analyzed  relative  to 
the  source,  is  shown  in  Figure  7b.  This  figure  show's  how  the  gust  front  can 
act  to  cut  off  the  low  level  inflow  to  the  updraft  if  it  propagates  well  ahead  of 
the  storm. 

A  simulation  was  also  run  with  a  source  speed  of  20  m/sec  and  a  shear  of 
4  X  10  sec  .  The  results  of  this  case,  MD7,  are  given  in  Table  1.  Since  the 
source  was  moving  very  fast  in  this  case,  the  density  current  did  not  begin 
propagating  away  from  the  source  until  it  had  achieved  a  sufficient  temperature 
difference  and  depth  for  its  density  current  speed  to  be  greater  than  the  source 
speed.  Its  close  proximity  to  the  source  for  a  considerable  length  of  time  re¬ 
sulted  in  substantially  lower  temperatures  in  the  head. 


The  density  currents  in  the  simulations  with  shear  are  moving  through  an 
environment  with  an  effective  ambient  wind  equal  to  the  average  winds  over  the 
depth  of  the  current.  Therefore,  Eq.  (4)  must  be  used  to  calculate  the  appro¬ 
priate  value  of  k.  For  cases  MDfi  and  MD7,  and  for  case  MD5,  which  had  a 
stationary  source  but  constant  ambient  wind,  it  was  found  that  the  best  agree¬ 
ment  was  obtained  when  0.  85  U  was  used  in  Eq.  (4)  rather  than  0.  82  U.  There¬ 
fore,  the  k  values  shown  in  Table  1  use  this  form.  The  other  values  of  k  shown 
in  Table  1  were  calculated  using  Eq.  (3).  The  average  k  value  for  the  M-type 
simulations  was  k  =  0.70. 

1.4  Moist  Simulations  With  M-tvpe  Source 

To  study  the  effect  arc  cloud  generation  has  on  gust  front  motion,  two  simu¬ 
lations  were  carried  out  in  a  moist  atmosphere.  The  environment  had  a  relative 
humidity  of  50  percent  in  the  lowest  1.  5  km  and  20  percent  above  3  km  with  a 
linear  transition  layer  between. 

The  results  of  these  two  simulations,  cases  MM1  and  MM2,  are  also  shown 
in  Table  1.  The  high  diffusion  required  for  the  M-type  source  appeared  to  sig¬ 
nificantly  affect  the  arc  cloud  liquid-water  distribution.  The  Q-type  source 
allows  a  much  smaller  diffusion,  which  is  less  likely  to  adversely  affect  the 
liquid  water  distribution.  Therefore,  although  the  results  shown  in  Table  1 
are  consistent  with  the  Q-type  simulations,  the  effect  of  the  arc  cloud  on  the 
gust  front  motion  will  be  discussed  in  terms  of  the  Q-type  source  results. 

5.  Q  TYPE  SOURCE  RESULTS 
5. 1  Orv  Results 

The  results  of  the  Q-type  simulations  are  given  in  Table  3.  The  dry  simu¬ 
lations  were  carried  out  for  a  stationary  source  using  several  source  strengths 
and  stabilities.  The  dry  simulation  results  using  the  Q-type  source  were  very 
similar  to  those  of  the  M-type  source.  The  average  value  of  the  Froude  number 
was  k  =  0.  75,  and  is  somewhat  larger  than  the  M-type  value.  This  is  probably 
a  reflection  of  the  smaller  eddy  viscosity  used  in  the  Q-type  simulations.  Al¬ 
though  the  flux  of  cold  air  tended  to  pulsate  in  the  Q-type  simulations,  the  flow 
in  the  vicinity  of  the  density  current  head  became  fairly  steady  after  10  to  15 
minutes  of  integration. 

Figure  8  shows  the  potential  temperature  field  for  case  QD5  at  20  min.  The 
density  current  produced  in  this  simulation  is  very  similar  to  the  one  shown  in 
Figure  5,  which  was  produced  by  the  M-type  source.  The  smoother  nature  of 


Figure  8.  Potential  Temperature  Field  for  Case  QD5  at  20  min 

the  source  region,  despite  the  smaller  diffusion  used  in  the  Q-type  simulations, 
is  evident  in  a  comparison  of  these  two  figures.  The  height  of  the  head  was 
somewhat  less  in  the  Q-type  simulations  than  in  the  M-type  for  a  given  temper¬ 
ature  difference  and  stability.  The  reason  for  this  result  is  not  clear. 

The  lifting  experienced  during  gust  frontal  passage  was  determined  by  com¬ 
puting  the  trajectories  of  parcels  originally  in  the  calm  air  ahead  of  the  front 
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using  the  method  of  Schlesinger.  The  lifting  profile  for  this  case  is  shown  in 
Figure  9.  It  was  found  that  even  though  the  density  current  is  a  shallow  phenom¬ 
enon,  lifting  is  experienced  by  parcels  throughout  the  depth  of  the  troposphere 

with  the  amount  of  lifting  decreasing  linearly  with  height.  This  result  is  of 
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direct  relevance  to  the  work  of  Matthews  who  was  forced  to  compare  several 
lifting  profiles  in  his  1-D  model  simulations  because  of  lack  of  guidance  in  this 
aspect  of  gust  front  passage. 

Given  the  stability  of  the  environment,  9e/9z  =  1.46°C/km,  the  lifting  of 
air  by  the  head  of  the  current  results  in  a  cooling  of  the  column  above  the  head, 
which  contributes  to  the  surface  hydrostatic  pressure.  Therefore,  Eq.  (2)  gives 
a  different  result  than  Eq.  (3)  since  it  considers  only  the  temperature  difference 
between  the  low  level  environment  and  the  current  itself.  For  a  nearly  neutral 
environment,  as  in  case  QD1  where  9()/9z  =  0.  2°C/km,  the  lifting  does  not 
result  in  much  cooling,  thus  the  values  of  k  computed  using  Eqs.  (2)  and  (3)  are 
nearly  the  same.  Table  3  lists  the  values  of  k  computed  using  both  Eqs.  (2)  and 
(3).  It  is  clear  that  the  simulated  density  currents  satisfy  the  density  current 
speed  equation  very  well  when  the  pressure  form  Eq.  (3)  is  used.  The  average 
value  of  k  computed  using  Eq.  (3)  is  also  very  nearly  the  same  as  that  given  by 
observations  of  gust  fronts. 

18.  Schlesinger,  R.E.  (1973)  A  numerical  model  of  deep  moist  convection: 

Part  A  prototype  experiment  and  variations  upon  it,  J.  Atmos.  Sci. 
30:1374-1391. 

19.  Matthews,  D.A.  (1981)  Observations  of  a  cloud  arc  triggered  by  thunder¬ 

storm  outflow,  Mon.  Wea.  Rev.  109:2 140-2 157 . 
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Figure  9.  Lifting  Experienced  by 
Parcels  During  Gust  Front  Pas¬ 
sage  for  Case  QD5 


5.2  Moist  Results 

The  results  for  simulations  in  a  moist  atmosphere,  with  low  level  humidities 
varying  from  50  percent  to  7  0  percent,  are  given  in  Table  3.  The  potential  tem¬ 
perature  distribution  for  one  of  these  simulations,  case  QM2  at  25  min  is  shown 
in  Figure  10a.  Figure  10b  shows  the  arc  cloud  location  in  the  liquid  water  dis¬ 
tribution  superimposed  on  the  potential  temperature  field.  For  clarity,  only 
the  center  20  km  of  the  40-km  domain  is  shown.  The  potential  temperature  dis¬ 
tribution  does  not  show  the  typical  head  structure  associated  with  density  cur¬ 
rents.  The  existence  of  the  head  is  evident,  however,  in  the  raised  isentropes 
at  upper  levels.  The  cloud  is  suppressed  behind  the  head  by  the  subsidence 
existing  there,  but  not  completely  dissipated.  All  of  the  moist  simulations 
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Figure  10a.  Potential  Temperature  Distribution  for  Case  QM2  at  25  min 


Figure  10b.  Liquid  Water  Field  for  Case  QM2  at  25  min  Superimposed  on  the 
Potential  Temperature  Field.  Liquid  water  contour  interval  is  0.  5  g/kg. 

Only  the  center  20  km  of  the  40  km  domain  is  shown 


exhibited  are  cloud  development  except  QM1,  but  only  QM4  had  precipitation 
reach  the  ground  under  the  arc  cloud. 

The  arc  cloud  is  produced  by  the  low  level  air  being  lifted  by  the  current 
above  its  condensation  level.  The  release  of  latent  heat  associated  with  this 
condensation  produces  a  warm  pocket  at  the  top  of  the  head,  which  results  in 
the  isentropic  structure  shown  in  Figure  10b.  The  air  over  the  head  has,  in 
fact,  been  raised  above  its  level  of  free  convection  and  is  warmer  than  the  envi 
ronment.  However,  deep  convection  does  not  develop  even  in  this  conditionally 


very  unstable  atmosphere.  It  is  felt  that  the  convection  is  suppressed  for  two 
reasons.  First,  the  air  lifted  by  the  current  is  quickly  raised  above  the  head 
and  then,  as  the  head  passes,  it  subsides  in  the  wake  region.  An  individual  par¬ 
cel  is  only  above  its  level  of  free  convection  for  a  short  period  of  time  before 
being  forced  downward.  Secondly,  the  frontal  boundary  above  the  head  is  a 
region  of  strong  shear,  which  will  inhibit  convection.  The  lack  of  deep  convec¬ 
tion  initiation  is  an  interesting  result.  Satellite  observations  of  gust  fronts  have 
shown  them  to  propagate  large  distances  producing  shallow  cumulus,  but  initiate 

deep  convection  only  after  they  intersect  with  another  outflow  boundary  that  pro- 
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vides  additional  lifting. 

The  production  of  the  warm  pocket  above  the  head  has  a  pronounced  effect 
on  the  speed  of  the  gust  front.  This  warm  air  decreases  the  surface  hydrostatic 
pressure  under  the  head  (although  the  weight  of  the  liquid  water  provides  a  slight 
increase)  and  therefore  reduces  the  pressure  gradient  force  driving  the  current. 
Thus,  the  moist  current  propagates  at  a  slower  speed  than  would  a  dry  current 
of  the  same  depth  and  temperature  difference.  This  fact  is  evident  in  Table  3 
through  the  comparison  of  the  dry  and  moist  simulation  results.  For  instance, 
comparing  QD4  and  QM2,  we  can  see  that  although  both  are  the  same  depth  and 
QM2  has  a  larger  temperature  difference  than  QD4,  QM2  propagates  slower  than 
QD4.  This  is  in  direct  conflict  with  Eq.  (2),  but  is  easily  explained  through 
Eq.  (3).  The  warm  pocket  above  the  head  in  QM2  results  in  a  smaller  surface 
pressure  difference,  and  therefore  a  slower  speed  of  propagation. 

Table  3  also  lists  the  calculated  updraft  velocity  at  the  front.  These  values 
compare  very  well  with  those  obtained  by  Mitchell  and  Hovermale.  11  This  up¬ 
ward  motion  is  primarily  accomplished  by  mechanical  lifting  even  in  the  moist 
cases.  Because  of  this,  for  a  given  stability,  the  updraft  speed  is  determined 
almost  entirely  by  the  speed  and  depth  of  the  density  current.  The  upward  mo¬ 
tion  is  suppressed  in  the  more  stable  environments,  as  would  be  expected.  If 
the  release  of  latent  heat  is  sufficient,  the  updraft  in  the  moist  cases  will  be 
augmented  by  the  buoyancy  of  the  parcels.  This  is  evident  in  a  comparison  of 
cases  QM2  and  QM4. 

6.  CONCLUSIONS 

This  report  has  presented  a  somewhat  different  version  of  the  gust  front 
speed  equation  based  on  the  surface  pressure  difference  between  the  density  cur¬ 
rent  head  and  the  environment.  It  was  shown  through  the  use  of  numerical  simu¬ 
lations  that  this  form  of  the  equation  gives  much  more  consistent  results  than  the 
other  commonly  used  forms.  It  is  also  easier  to  apply  to  observations  since  only 
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the  surface  values  of  the  thermodynamic  variables  are  needed.  Twenty  observed 
gust  front  cases  were  re-analyzed  in  terms  of  this  equation  and  the  appropriate 
internal  Froude  number  was  found  to  be  k  =  0.79. 

This  study  also  showed  that  accurate  representation  of  the  thunderstorm  out¬ 
flow  and  gust  front  can  occur  in  a  model  with  relatively  coarse  resolution  in  the 
vertical.  This  is  an  important  consideration  in  numerical  thunderstorm  models, 
which  have  a  similar  resolution  but  depend  on  gust  front  effects  for  maintenance 
of  the  storm.  It  should,  however,  be  noted  that  proper  treatment  of  the  surface 
drag  was  found  to  be  essential,  especially  at  this  relatively  coarse  resolution. 

Simulations  with  a  moving  cold  air  source  shed  light  on  the  gust-front  de¬ 
coupling  process.  Once  the  gust  front  has  propagated  a  sufficient  distance  from 
its  source  region,  its  speed  is  determined  by  its  local  characteristics.  How¬ 
ever,  the  source  region  can  have  a  great  deal  of  influence  on  the  outflow  if  the 
ambient  wind  and  source  speed  act  to  keep  the  gust  front  close  to  the  source. 

This  indicates  that  modulations  in  the  thunderstorm  strength  and  speed  may  act 
to  retain  the  gust  front  in  a  position  which  is  beneficial  to  the  storm.  Any  pos¬ 
sible  feedback  mechanism  that  could  control  the  modulations  must  keep  the  gust 
front  from  propagating  beyond  its  decoupling  distance  in  order  to  produce  a  long- 
lived  storm. 

The  warm  pocket  produced  above  the  head  during  arc  cloud  formation  tends 
to  result  in  a  slower  gust  front  than  a  corresponding  outflow,  which  did  not  pro¬ 
duce  an  arc  cloud.  This  decrease  in  speed  is  a  result  of  the  decreased  surface 
hydrostatic  pressure  rise.  Deep  convection  is  not  initiated  by  the  gust  front, 
even  when  parcels  are  lifted  above  their  level  of  free  convection  during  gust 
front  passage,  indicating  that  an  additional  trigger  mechanism  may  be  necessary. 

This  study  has  produced  some  very  interesting  results  concerning  the  effects 
of  ambient  winds  and  moisture  on  density  current  motion.  Many  of  these  results 
have  direct  bearing  on  the  speed,  depth,  and  intensity  of  thunderstorm  outflows 
and  their  associated  gust  fronts.  Further  work  is  required,  however,  to  truly 
quantify  these  effects.  The  result  that  the  lifting  is  experienced  throughout  the 
depth  of  the  troposphere  indicates  that  results  from  an  anelastic  model  should 
be  included  in  this  study.  Further,  a  more  complete  series  of  moving  source 
and  sheared  environment  studies  should  be  carried  out  in  a  more  stable  envi¬ 
ronment,  where  the  density  current  would  have  a  more  realistic  depth.  This 
shallower  depth  would  require  a  higher  resolution  in  the  vertical  (say,  200  m) 
in  order  to  allow  the  proper  interaction  of  the  current  with  the  shear.  The  ef¬ 
fects  of  the  shear  may  be  considerably  reduced  for  shallow  currents. 
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